Continuum model for polymers with finite thickness 
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Abstract 

We consider the continuum limit of a recently-introduced model for dis- 
cretized thick polymers, or tubes. We address both analytically and numeri- 
cally how the polymer thickness influences the decay of tangent-tangent cor- 
relations and find how the persistence length scales with the thickness and 
the torsional rigidity of the tube centerline. At variance with the worm-like 
chain model, the phase diagram that we obtain for a continuous tube is richer; 
in particular, for a given polymer thickness there exists a threshold value for 
the centerline torsional rigidity separating a simple exponential decay of the 
tangent-tangent correlation from an oscillatory one. 

Experimental studies of biopolymers have constantly stimulated the search for 
schematic models apt for reproducing the observed kinetic and thermodynamic be- 
haviour. In recent years two types of biopolymers have attracted most of these efforts: 
DNA and proteins. The interest in the former has been sparked by the introduction of 
single-molecule experiments which probed the elastic response of DNA upon stretch- 
ing. 1 Considerable progress in the rationalization of these experiments was made 
thanks to the worm- like- chain (WLC) model 2,3 where the polymer is described as a 
continuous centerline possessing an effective bending 2 and/or twisting 4 ' 5 rigidity. For 
protein modelling, instead, one of the goals is to capture the main physico-chemical 
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forces responsible for driving the folding process towards the native state. Typical 
coarse-grained models adopt sophisticated energy functional (often with hundreds 
of parameters) in order to reproduce the observed variety of protein folds. It has 
been recently argued, however, that the overburdening of the energy function can be 
avoided by modelling explicitely the intrinsic thickness of proteins. 6-9 As we discuss 
later this is achieved through the introduction of suitable three-body interactions 
among triplets of points constituting the polymer centerline. 8,10-13 It is physically 
appealing that the thick-chain model has proved valuable also for the case of DNA 
in applications ranging from the characterization of knotted DNA molecules 14 to the 
thermodynamics of DNA packaging inside viral capsids and DNA elastic response. 15 
By necessity, all numerical implementations of these models rely in the discretiza- 
tion of the polymer centerline into a succession of beads whose "natural" spacing is 
typically suggested by the intrinsic granularity of the polymer itself (e.g. the sepa- 
ration of consecutive C a 's for proteins or the base-pair spacing in dsDNA). From a 
theoretical perspective it is therefore desirable to characterize the thick chain model 
in the continuum limit, where the bead spacing tends to zero (analogously to the 
WLC limit of the Kratky-Porod model). This continuum limit has, so far, been 
considered only for characterizing the limited repertoire of ideal knots. Motivated 
by the potentially-wide range of applicability of the thick-polymer model, in this 
Letter we take the perspective of addressing the statistical mechanics of general- 
ized thick-polymer models in the continuum limit. In particular, we introduce in 
the Hamiltonian a penalty for the geometrical torsion of the tube centerline and, 
initially, consider the constraints induced by the finite polymer thickness only at a 
local level, a simplification usually adopted to allow analytical progress. 4 ' 5 From the 
exact analysis it emerges that accounting for the centerline torsional rigidity term 
(1) allows to get a finite persistence length in the limit of a continuous thick polymer 
and (2) introduces a novel feature in the behavior of the tangent-tangent correlation 
function, namely the presence of a threshold value for the torsional rigidity which 
separates a monotonic decay from an oscillatory one. Finally, numerical Monte Carlo 
simulations are employed to show that this behavior persists also when the tube con- 
straint is enforced at non local level. These findings highlight the rich behaviour of 
models where the thickness is treated explicitly. As a comparison we consider the 
case of a WLC in the presence of penalty for the centerline geometrical torsion. It 
is found that this model exhibits either a simple exponential decay or an oscillatory 
one (with singular behaviour in the continuum limit) depending, respectively, on 
the absence or presence of the torsional rigidity but independently of its strength. 
In relation to the behaviour observed here, it should also be noted that Panyukov 
and Rabin 16 have considered, in place of the whole equilibrium ensemble, a rod-like 
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Figure 1: Frenet reference frame for a discrete bead model. Note that 0i e [0,7r] is 
defined for 1 < i < N — 2, whereas 0, G [0, 2tt] is defined only for 2 < i < N — 2. 



chain fluctuating around a stress-free helical conformation. Under these conditions 
they could observe a change from an oscillatory to a simple exponential decay by 
increasing the fluctuation strength. 

We model a polymer chain by means of a set of N consecutive beads, {ro---rzv-i} 
connected by bonds of fixed length, a. The succession of beads constitutes the 
centerline for our thick polymer. We shall denote with A and n t the thickness of the 
chain and the torsional rigidity, respectively. Although we shall first focus on the 
case K t = we will develop a formalism general enough to be used also in the WLC 
with torsional penalty. By analogy with the Frenet reference frame for continuous 
curves, 17 we define an orthonormal set associated to each bead, formed by the local 
tangent, ti = (r i+ i — r^)ja, binomial, 6, = ti A tj_i/|ti A tj_i| and normal, hi = ti A b{ 
(see Fig. ^i). It is possible to write recursion equations relating the reference axes 
for bead i in terms of those for bead i — using the polar angles 9i and <pi as in Fig. 
Eh: 

h = cos0i&i_i - sin0ini_i , 
< ti = sin 9i sin + cos 9^^ + sin $i cos 4>ihi^i , 

hi = cos 9i sin 4>S%~\ — sin + cos f cos 4>ihi_i . 

Quite generally the joint probability distribution of angles, V{9\, 6*2, 02, #3, 03, • • ■), 
resulting from the canonical average, will depend on the whole ensemble of inter- 
actions including the steric ones. However, for the simplified case where the effects 
of the polymer thickness are treated only locally (as for twist and bending rigidity) 
then V factorises in terms of the probability distributions p(6 l j,0j) for each pair of 
angles 9 h 4>i- 

N-2 

V(9 1 , 02, 2 i 03, 03, ■■■) = Pe{0i) II ■ (!) 

i=2 

Since we are considering a uniformly-thick homopolymer, the same probability distri- 
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bution, p (9, 0), is involved for all beads. In the following the averages weighted with 
p will be denoted as (•), while those calculated with respect to V will be written as 
(-)-p- The factorization leads to a straightforward characterization of the decay along 
the chain of expectation values such as j\ = (bi ■ x)p, g,i = (t« • x)p, hi = (hi ■ x)-p, 
where x could be, e.g., t±, b\, hi. In fact, /, g and h at location i + 1 are obtained 
from those at site i by the application of the following transfer matrix: 

/ (cos0) — (sin0) \ 

T= (sin6»sin0) (cos0) (sin6>cos0) . (2) 
\ (cos 9 sin 0) — (sin#) (cos 9 cos 0) J 

If the eigenvalues of T are real, the decay of /, g and h will be monotonic, while if 
they are imaginary there will be an oscillatory modulation. 

We now consider two further simplifying assumptions: (i) the 'bond' angle 9 
and the 'dihedral' angle contribute independently to the probability distribution 
p{9,4>) = pe(6)p ( f > ((f)), and (ii) the system is invariant for chirality flipping p«j!,(0) = 
p^(—(f)). In this case, (sin0) = 0, and the transfer matrix T becomes block diagonal 
with an eigenvalue Ai = (cos0), so that \b n -bi) = A™ -1 = exp [— a (n — 1) /£, 
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decays exponentially with a correlation length = — a/ln(cos0). The remaining 
two eigenvalues are the solutions of the second order equation A 2 — b\ + c = 
with b = (cos 9} [1 + (cos0)], c = (cos0) [(cos#) 2 + (sin#) 2 ]. The relevant quantity 
which discriminates between different decay properties is T = b 2 — 4c . If T > 
0, the two solutions A2,3 = (±VT + b)/2 are real and the correlation function for 
tangent /normal vectors decays exponentially to zero, with the correlation length 
£ t = —a/ In A 2 being controlled by the largest eigenvalue. If T < 0, the two solutions 
are complex conjugate, and the tangent-tangent correlation function exhibits an 
oscillatory decay: 

(i n+1 -i 1 ) = C ° S[a/X0 + an/x] e— /*, (3) 
cosxo 



where £ t = —2a/ lnc, x = a l arctan {\/— r/&J , and xo depends on initial conditions. 

Within this general framework we now consider different specific examples. We 
first focus on the case of a WLC subject to a penalty for the geometrical torsion; the 
corresponding Hamiltonian is: 

where Kb, K t , are the bending and torsional rigidity, respectively, defined in such a 
way to get back the usual Hamiltonian in the continuous limit, a — > 0. We emphasize 
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the fact that the Hamiltonian (JH) is entirely specified by the centerline geometry The 
approach therefore differs in spirit from those used to model the elasticity of rod-like 
chains. In these contexts, starting from a stress-free refecence configuration, one 
introduces a "material reference frame" which is used to keep track not only of the 
deformations of a given reference frame, but also of the twist around the centerline. 17 
Although this latter information is clearly not available in the model of eqn. (J3| it 
is worth to considering energy-functions relying uniquely on the knowledge of the 
centerline. In fact, for the important class of biopolymers constituted by proteins, it 
is well known that the knowledge of a protein's centerline (the C a trace) and sequence 
composition allows to reconstruct the whole protein structure with high accuracy. 

For the case of Hamiltonian (0J) the probability distributions for bond and dihedral 
angles taking, of course, into account the inverse temperature (3 is found to be: 
Po (9) = sin^exp (3^ cos 9 ; p^ (0) = exp (3^ cos0 . All averages appearing in the 
transfer matrix elements of eq. (J2J) can formally be expressed by means of modified 
Bessel functions which, in turn allow to identify the boundary, V — 0, separating 
the oscillatory from the monotonic decay of tangent correlations, see Fig. In 
the continuum limit, a — > 0, the angles 9 and <fi contributing significantly to the 



average come from a region centred around zero and of width JajftKi, and Ja//3rz t , 
respectively and the equation for the boundary is a/ /3k& ~ (a/ (3n t ) 2 /8ir. This implies 
that, for any finite value of the torsional rigidity, in the continuum limit the tangent- 
tangent correlation function always decays in an oscillatory way. The period of the 
oscillation is proportional to yfa K h while the decay length, £ t , is independent of a, 
^t 1 = l/(&(3K t ) + (1 — 7r /4)/(/?K fe ). It is therefore apparent, that in the continuum 
limit, a — > 0, the oscillation period becomes smaller and smaller, denoting a singular 
behaviour of the chain. This is reminescent of the singular behaviour of rod-like 
chains which, in the continuum limit, exhibit plectonemes on smaller and smaller 
scales. 4 ' 5 From Fig. |21 it is apparent that only if «t is exactly zero, one remains in the 
vanishingly small region of monotonic exponential decay when the continuum limit 
is taken. In this case one recovers the WLC case with persistence length £ t ~ (3^. 

We consider now the case of a polymer chain describing a thick self-avoiding tube 
of uniform cross-section. The finite thickness A of the polymer impacts on two dis- 
tinct conformational features. First, it constrains the local radius of curvature to 
be not less than A. 10 ' 11 In addition, there is also a non-local effect since any two 
portions of the tube, at a finite arclength separation, cannot interpenetrate. 10 ' 11 In 
traditional beads-and-strings models it is only this second effect that is taken into 
account through a pairwise potential with a hard-core repulsion. Interestingly, one 
needs to go beyond pairwise interactions to account for the above mentioned effects 
in discretized polymer chains. 7, 8 ' 12 ' 13 In fact, the requirement on the local radius 



5 



0.1 




0.08 



^ 0.06 
^ 0.04 
0.02 











0.25 



0.5 



0.75 



a/pk, 



Figure 2: Boundary separating the oscillatory from the exponential decay of the 
tangent-tangent correlation. 

of curvature can be enforced by finding the radii of the circles going through any 
consecutive triplet of points and ensuring that each of them is greater than A. The 
non-local effect can be addressed within the same framework by considering the min- 
imum radius among circles going through any non- consecutive triplet of points is also 
greater than A. In summary the finite thickness, A, of the discretized tube requires 
that the radii of the circles going through any triplet of distinct points have to be 
greater than A (see Fig. |SJ). 11-13 In the present context, we are interested mainly in 
the local thickness effects. Therefore, we will consider the following reduced Hamilto- 
nian in the absence of torsional rigidity, involving only local three-body constraints: 
= J2i V{Ri-i,i,i+i)i where R^ is the radius of the circle going through the beads 
i, j, k, and the potential V(R) is oo, if R < A, and otherwise. From simple 
geometric considerations, the local tube constraint can be expressed in terms of the 
bond angle distribution by imposing pe (9i) — if 9\ > 2arcsin (^H- One can see 




Figure 3: Sketch of a curve which is a viable centerline for a tube of thickness A. 
The radii of the circles through any triplet of points, are not smaller than A. 

that the thickness A plays a role similar to the bending rigidity Kb, in that they both 
induce the chain to be locally more straight. Yet, the scaling behavior in the a — > 
limit is different in the tube case, since the range of 9 angles most contributing to 
the average has now width ^ (instead of ^a/Pn b ), yielding (cos 9) = 1 — a 2 /4A 2 and 
(sin#) ~ 2a/3A. Since T = (cos9) 2 ~ 1 — a 2 /2A 2 , the tangent-tangent correlation 
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function decays exponentially. This is similar to the WLC case in the absence of 
torsional rigidity, but in the tube model the correlation length scales in a different 
way £f ~ — . This is a pathological behavior in the continuous limit, since the 
correlation length diverges as a — > 0. The recipe by which the tube constraint is im- 
plemented for a discrete chain ends up in preferentially sampling straight continuous 
lines. It is natural to associate such ill behavior to the degeneracy in the choice of 
the Frenet frame which arises for a straight line conformation. 18 This can in fact 
be cured by adding a torsional rigidity term to the tube constraint, constraining the 
unphysical fluctuations of the binormal vectors around straight line conformations. 
The Hamiltonian for this rod-like thick polymer is: 



(5) 



In this case, (cos< 
whereas we get Y • 



~ 1 — a/2/3 /% in the a — > limit, which implies = 2(3 n t 
36/3 " K ^ A2 [9A 2 — 64/3 2 /t 2 ]. Thus, there are two different regimes, 



oscillatory decay if n t > k% = ~~, and monotonic decay if K t < k*[. In the latter 
case the persistence length associated with the tangent-tangent correlation function 
is controlled by both the tube thickness and the torsional rigidity: 



9A 2 



1 + 



\ 



'8/W 
3A 



(6) 



In the former case, the correlation length depends instead only on the torsional 
rigidity, whereas the oscillation period depends also on the tube thickness: 



it = , x = — 



3A N 



-1/2 



(7) 



It can be seen from eq. that by increasing the torsional rigidity the tangent- 
tangent correlation length £ t , initially controlled by the thickness A decreases until 
the threshold k* is reached. Above such threshold torsional rigidity takes over and 
the tangent-tangent correlation length becomes equal to twice the binormal-binormal 
correlation length, but the thickness signature remains in the oscillatory behavior of 
the tangent-tangent correlation function and in the period x- The previous analysis, 
is in good semi- quantitative agreement with data from MC simulations on the full 
model (i.e. where non-local effects are retained), as visible in Fig. 0] The simulations 
were performed on chains of 128 equispaced beads through the Metropolis acceptance 
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Figure 4: Tangent-tangent correlation, C(s), as a function of the arclength separa- 
tion, s. The dashed and solid lines refer respectively to the analytical results from 
eq. (J2J) and to the MC simulations (maximum dispersion pa 0.01). 

of crankshaft and pivot moves. The tangent-tangent correlations were measured by 
sampling structures at intervals greater than the autocorrelation time. 

To conclude, we have shown how the recently introduced thick-polymer model 
can be regularized in the continuum limit by introducing a term penalizing the geo- 
metrical torsion of the centerline. For any given value of the polymer thickness there 
exists a torsional-rigidity threshold separating the monotonic decay of the tangent- 
tangent correlation from the oscillatory one. This highlights the rich behaviour of 
thick-polymer models which combine features previuosly observed in distinct poly- 
mer models, such as the worm- or rod-like chains. The wide use of the latter in the 
context of single-molecule experiments opens the possibility to use the physically- 
appealing perspective of semi-flexible thick polymers to interpret the behaviour of 
biopolymers. 

We thank J. R. Banavar, T. X. Hoang and F. Seno for useful discussions and 
acknowledge support from COFIN MURST 2003, FISR 2001, INFM and EPSRC. 
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